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Motivation 


► Computational Aeroacoustics (CAA) codes are increasingly 
used to simulate complex physics, make design decisions, etc. 

► All software has bugs[l]. How can we know we’re getting the 
right answer out of our CAA code? 
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Code Verification 


► Code verification: “are we solving the equations right?” All 
about the math/numerics. Goal is to detect bugs. 

► Gold-standard of code verification: the order-of-accu racy test. 
Does the error behave as we expect, i.e., converge at the rate 
of the order-of-accuracy of the code’s schemes, p? 

► e « Ah p 

► So, we need the code’s error, i.e., the difference between the 
code’s solution and a reference solution. 

► Where do we get a reference solution? 

► Method of Exact Solutions (MES) 

► Method of Manufactured Solutions (MMS) 
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Reference Solutions for Code Verification: State of the Art 

Method of Exact Solutions (MES) 

Look in a textbook (or paper, etc.) for a solution to the PDE the 
code solves. 

► Advantages: Simple, no modifications to PDE code. 

► Disadvantages: available solutions tend to be either tricky to 
evaluate numerically (infinite series, integrals, etc.), or too 
simple (terms go to zero, etc.), or both. 

Method of Manufactured Solutions[2, 3] (MMS) 

Start with the solution you want, then change the PDE to make it 
work. 

► Advantages: General and flexible — user chooses the solution. 

► Disadvantages: Source terms must be added to the PDE, and 
thus the code, which can be quite complicated for non-trivial 
PDEs. 


The Ideal Code Verification Method 


► The complication of MMS source terms appears to limit the 
popularity of MMS. 

► Characteristics of an ideal code verification method? 

► Flexibility of MMS: user has control over the form of the 
solution. 

► Unobtrusiveness of MES: no modifications of PDE, or code. 

► External Verification Analysis (EVA) is intended to fulfill these 
requirements. 

Goals of this work 

► Show how EVA can provide a reference solution for the 
nonlinear Navier-Stokes equations suitable for code 
verification, 

► Use EVA to verify a high-order, viscous Computational 
Aeroacoustics code. 


EVA: An Approximate Solution to the Cauchy Problem 

► Express PDE as a Cauchy (initial value) problem, 
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► Use a Taylor series to approximate u(x, t) as 
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► The EVA approach: use the above truncated series combined 
with an analytic initial condition as a reference solution for 
code verification. 

► Cauchy-Kowalewski (CK) recursion: repeatedly differentiate 
governing equation, eventually expressing. . . 

d n u d m u 
dt n ^ dx m 

► Differentiate the analytic initial condition exactly. 

► Details in the paper! 


A Recurrence Relation for (Derivatives of) the 
Navier-Stokes Equations 

► Recurrence relation for the continuity equation: 

Qb+d+n+l^ = /b\ y / d \ " /n\ 

dx b dy d dt n+1 ^ \a) ^ \c ) \m J 
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► Using the specific volume a = - (makes the math easier). 

► Not as bad as it looks. Like Dyson and Goodrich[4], using 
Leibniz rule to differentiate products. 


The Initial Condition 


► Final piece of the puzzle: initial condition (1C). EVA uses an 
analytic 1C. 

► Since EVA solves a Cauchy (read: initial value) problem, no 
particular boundary conditions (BCs) are enforced. So an 1C 
with a “footprint" helps us avoid problems with the PDE 
code’s BCs. 

► The 1C currently implemented in EVA tool: 
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► Need to be able to calculate arbitrary-order derivatives of the 
initial condition. 

► Benefit from analytic derivatives: EVA solution is completely 
unaffected by grid/mesh used! 


The Computational Aeroacoustics Code: BASS 

► Broadband Aeroacoustic Stator Simulator from NASA Glenn 

► High-order, parallel, block-structured finite-difference 
Computational Aeroacoustics solver with explicit, optimized 
time-marching 

► Inviscid or Viscous, 2D or 3D 

► Spatial differencing schemes 

► E_2: Explicit 2 nd -order 

► E_6: Explicit 6 th -order 

► DRP: Tam & Webb Dispersion Relation Preserving[5] 

► C_6: Hixon's[6] prefactored form of Lele’s[7] compact 6 th -order 

► Time marching schemes 

► RK4L, RK5L: Jameson[8] four- and five-stage Runge-Kutta, 

► RK56: Stanescu and Habashi[9] five/six stage Runge-Kutta, 

► RK7S, RK67: Allampalli et a I . [10] High-Accuracy Large-step 
Explicit Runge-Kutta (HALE-RK) seven and six/seven-stage 
schemes. 

► Modified BASS to use constant viscous properties. 


Verification Process 


1. Identify (at least) formal order-of-accuracy of CAA code’s 
numerical schemes, 

2. Choose test case parameters, 

3. Run EVA, 

4. Run CAA code with range of “discretization measures" (grid 
spacings or time step sizes), 

5. Calculate error and some error norm, 

6. Calculate error norm convergence rate, compare to CAA 
code’s formal order-of-accuracy. 


Spatial Test Case Parameters 


► Skewed, rotated, and “perturbed” 3D grid series. 

► RK7S time marching scheme with a single very small time 
step At = 0.0005 (far-field acoustic velocity is 1). 

► Initial condition has an approximate wavenumber 

k = ^5 = 87r, which, when combined with the grid spacings 
chosen, give a range of points-per-wavelength from 2 to 50. 
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Test Case Parameters: Grid 


Test Case Parameters: Gaussian Centers 
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Test Case Parameters: p Initial Condition 
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Test Case Parameters: pw Initial Condition 
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BASS Spatial Verif.: p , fa, AAx p Assumption 
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BASS Spatial Verif.: Multiblock Grid for Parallel Runs 
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BASS Spatial Verif. : p , fa, AAx p Assumption, Parallel 




BASS Spatial Verif.: pv, C_6 Serial vs. Parallel 



Temporal Test Case Parameters 


► Grid and initial condition identical to the spatial test cases. 

► Used E_6 scheme for spatial differencing. 

► Marched to a final time level of t — 0.05 with an EVA target 
truncation error of 10 12 , which required a 26 th -order Taylor 
series. 

► Combination of time step size and initial condition gives a 
steps-per-period range of 5 to 200. 
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BASS Temporal Verif. : p, fa, AAt p Assumption 
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BASS Temporal Verif.: p, fa, AAt p + B Assumption 
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BASS Temporal Verif. : RK4L Low vs. High Amp. 
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Conclusions 


► EVA provides a flexible method for obtaining a reference 
solution suitable for code verification 

► Extended to support 3D nonlinear Navier-Stokes. 

► No modification to CAA code required. 

► EVA doesn't care about the quality of the grid (just like an 
exact or manufactured solution). 

► EVA was used to verify the spatial and temporal schemes in 
BASS, a high-order CAA solver. 

► Each scheme eventually converged at the expected rate, 
strongly verifying the code. 

► Overall, error and convergence rates corresponded well to 
linear analysis of the schemes (see paper). 

► Needed more general AAt p + B assumption when calculating 
convergence rate for the high-order time-marching schemes. 

► Saw linear/nonlinear behavior for RK4L scheme (4 th -order 
linear, 2 nd -order nonlinear). 
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Future Work 


► Current limitations of the EVA approach to code verification: 

► Boundary conditions, 

► Non-constant properties (i.e., viscosity and thermal 
conductivity), 

► Finding the asymptotic range. 

► More codes(?) 
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The End 


Thanks! 

This work was supported by the NASA Advanced Air Transport 
Technology (AATT) Project. The authors would like to thank 
Dr. Edmane Envia of the NASA Glenn Research Center, who was 
the technical monitor for this work. 
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MES Example 

From Knupp & Salari [1 1] , the exact solution to unsteady 
homogeneous heat conduction in a solid sphere is 


-j- oo oo n ^ — cxAppt 
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Takeaway 

Exact solutions are simultaneously too complex, too simple. 
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What is the Method of Manufactured Solutions[2, 3]? I 


► The “backwards approach:" start with the solution you want, 
then change the PDE to make it work. 

► Simple example with the linear advection equation: 

1. The manufactured solution: 


u(x, t) = sin(x — bt) 


2 . 

3 . 


The PDE: 


du du 


Get the MMS source term by putting the manufactured 
solution in the PDE: 


S(x, t) 


du du 
dt 3 dx 

[sin(x - bt)] t + a [sin(x - bt)] x 
( a — b) cos(x — bt). 
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What is the Method of Manufactured Solutions[2, 3]? II 

4. Add S(x, t ) to the PDE to make u an exact solution: 

du du . . . . , , . 

— + a— = S(x, t) = (a - fa) cos(x - fat). 

5. Use u as the reference solution. 

► Advantages 

► General: applicable to any PDE. 

► Flexible: user chooses the solution, so all aspects of PDE code 
can be tested 

► Disadvantages 

► MMS source term can be complicated (but computer algebra 
system software [CAS] helps) 

► PDE code must allow for user-specified distributed sources, or 
be made to. 
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MMS in the “Real World” 


► Surveyed a recent volume of the Journal of Computational 
Physics[12] and Computers & Fluids[13] (total of 69 articles) 
► Was code verification performed? If so, how? 



Cauchy-Kowalewski Recursion: What is it? 


► Cauchy-Kowalewski (CK) recursion: fancy math technique for 
getting the temporal derivatives needed for the Taylor series 
coefficients by repeatedly differentiating the original PDE. 

► Used in the proof of the Cauchy-Kowalewski theorem, a local 
existence and uniqueness theorem for PDEs. 

► Occasionally used to integrate PDEs numerically: 

► Dyson & Goodrich's MESA [14, 4, 15] 

► Castro & Toro's ADER [16, 17] 


CK Recursion: How it works 


► Goal: to express in terms of 

► Two “phases" CK recursion: 

1. Find by taking t-derivative 

2. Look at finding all unknown mixed-x-t derivatives by 
taking x-derivatives of known expressions, repeating until only 
pure x-derivatives remain. 

► See paper for all the (very exciting!) details. 
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CK Recursion: Viscous Burgers Dependency Graphs 
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Differentiating Governing Equations with the Leibniz Rule 


► Problem: CK recursion requires repeated differentation of the 
governing equation (here, 3D nonlinear Navier-Stokes). 

► The terms in the flow equations that give us trouble are 

products etc.). 

► Is there a pattern? 


► 
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Differentiating Governing Equations with the Leibniz Rule 


► The Leibniz rule: generalization of the product rule from 
calculus: 

k=0 ' ' 


► 

(”) : (n, k) binomial 

coefficient (Pascal’s 

triangle 

So, 
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0 f-g) i2) 
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► 

( f-g) {3) 

= fwg + 

3 f{2)g(l) +3 f(l) g {2) 

+ te (3) 


► Approach: As in Dyson[4], repeatedly use Leibniz rule to 
differentiate with respect to each variable (x, y, z, t ), giving 
recurrence relation for all the derivatives needed for CK 


recursion. 


Recurrence Relation for the Navier-Stokes 

► The x-momentum equation: 
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Order-of-Accuracy Calculation 


Two different approaches used here: 

► Simplest approach: assume 

e; = & - O « Ah? 

Two unknowns (A, p ), requires error from two CAA runs. 

► Slightly more complicated: assume 

e; « Ah? + B 

which has three unknowns (A, B, p) — needs three CAA runs. 
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Spatial Verification: Fourier Analysis 


► At a minimum, need to know scheme's order-of-accuracy to 
verify with EVA (or anything else) and the order-of-accuracy 
test 

► But a more detailed understanding of a scheme's performance 
helps with setting up test cases and interpreting results 

► Fourier analysis for finite differencing schemes: use scheme to 
differentiate a single Fourier component f(x ) = e lkx , then 
compare to the exact value. 


FD Schemes: Fourier Analysis Error and Conv. Rate 




Temporal Verification: Fourier Analysis 


► Fourier analysis for time-marching schemes: use the scheme 
to integrate the ODE 


du 

— = — ILdU 

dt 


which has the exact solution 


u(t) = u 0 e~ iujt 


where uq = u( 0). 


Temporal Schemes: Global Error and Convergence Rate 
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